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The PHENIX experiment has measured the suppression of semi-inclusive single high transverse 
momentum tt^'s in Au+Au collisions at ^sjvjv ~ 200 GeV. The present understanding of this 
suppression is in terms of energy-loss of the parent (fragmenting) parton in a dense color-charge 
medium. We have performed a quantitative comparison between various parton energy-loss models 
and our experimental data. The statistical point-to-point uncorrelated as well as correlated sys- 
tematic uncertainties are taken into account in the comparison. We detail this methodology and 
the resulting constraint on the model parameters, such as the initial color-charge density dN^ /dy, 
the medium transport coefficient (q), or the initial energy- loss parameter eo- We find that high 
transverse momentum tt" suppression in Au-f-Au collisions has sufficient precision to constrain these 
model dependent parameters at the ± 20-25% (one standard deviation) level. These constraints 
include only the experimental uncertainties, and further studies are needed to compute the corre- 
sponding theoretical uncertainties. 



PACS numbers: 25.75.Dw 

I. INTRODUCTION 

Heavy ion collisions at very high energy are of inter- 
est due to the formation of a novel partonic medium 
approximately the size of a large nucleus, but with an 
energy density exceeding that of normal nuclei by con- 
siderably more than an order of magnitude. At such high 
energy densities, it is believed that quarks and gluons are 
no longer confined in hadrons, but may be constituents 
of a quark-gluon plasma with characteristics of a near- 
perfect fluid (for a detailed review see Experiments 
at the Relativistic Heavy Ion Collider (RHIC) have al- 
ready demonstrated that a very hot and dense, strongly 
interacting medium is created in Au-f Au collisions at 
^/JWW^20Q GeV 0]. The goal is now to quantitatively 
determine the properties of this medium. 

Important properties of the medium include the den- 
sity of color-charges as well as the exchange of transverse 
momentum between parton probes and the medium. In 
rare events, in addition to the creation of the medium, 
there can also be a hard scattering (high-Q^ process) 
between the colliding partons that, at leading order, 
sends two high-energy quark or gluon partons in oppo- 
site transverse directions. These high-energy partons can 
be utilized to probe both the color-charge density of the 
medium and the coupling strength between the parton 
and the medium. There are various calculational frame- 
works for modeling these interactions (for a detailed re- 
view see 3). 

In this paper, we consider four specific calcula- 
tions of parton energy- loss (discussed below): the Par- 
ton Quenching Model (PQM) [1], the Gyulassy-Levai- 
Vitev (GLV) model the Wicks-Horowitz-Djordjevic- 
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Gyulassy (WHDG) model i, and the Zhang-Owens- 
Wang- Wang (ZOWW) model 0. We detail a quanti- 
tative method of assessing the sensitivity of the latest 
measurements to the input parameters of these models 
that characterize the initial parton density or medium 
transport coefficients. 



II. EXPERIMENTAL RESULTS 

During the 2004 data-taking period at the Relativistic 
Heavy Ion Collider, the PHENIX experiment recorded 
an integrated luminosity of 0.24 nb~^ in ^snn=200 GeV 
Au-f Au collisions, which extends the measurement of 
to much higher pT than previous data sets allowed. The 
results and further details of this measurement are given 
in [9] . A brief description is given below. 

The PHENIX experiment measures tt^'s via the two- 
photon decay mode with two types of highly segmented 
(A77 X A<j) « O.Ol X O.Ol) electromagnetic calorime- 
ters (EMCal) located at the radial distance of approxi- 
mately 5 m from the vertex Q . One is a lead scintillator 
sampling calorimeter (PbSc) , which covers the geometri- 
cal acceptance of I77I < 0.35 and A0 = 3/47r. The other 
is a lead glass Cerenkov calorimeter (PbGl), whose ge- 
ometrical coverage is |?7| < 0.35 and Acj) = 7r/4. The 
energy resolution of the PbSc and PbGl calorimeters as 
determined from test beam measurements are given by 
8.1%/^£;(GeV)©2.1% and 5.9%/^£;(GeV) ® 0.8%, re- 
spectively. The energy calibration of the EMCal mod- 
ules is based upon the measured position of the 7r° mass 
peak, the deposited energy of minimum ionizing particles 
which traverse the calorimeter, and the ratio of energy 
to momentum which is expected to be about 1 for elec- 
trons identified by the Ring-Imaging Cerenkov detector. 
The systematic uncertainty on the energy scale is ^ 1 %, 
which corresponds to ~ 7-12 % uncertainty on the invari- 
ant tt" yield over the px range of the measurement. 
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Neutral pions were reconstructed in their 7r° 77 
decay channel. Photon candidates are identified by ap- 
plying particle identification cuts based mainly on the 
shower shape. The invariant mass for all photon-pair 
combinations within one event that satisfy cuts on the 
energy asymmetry — E^2\/\E~fi + Ey2\ < 0.8 were 
calculated in bins of px- The combinatorial background 
is determined by combining into pairs uncorrelated pho- 
tons from different events with similar centrality, reaction 
plane, and vertex location. 
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FIG. 1: (Color online) The tt° nuclear suppression factor -Raa as a 
function of transverse momentum for 0-5% Au+Au collisions at 
^/slvlv— 200 GeV. Point-to-point uncorrelated statistical and 
systematic uncertainties are shown as uncertainty bars. Correlated 
systematic uncertainties are shown as gray boxes around the data 
points. The global scale factor systematic uncertainty is ±12%. 

The raw tt" yield was obtained by integrating the mass 
peak region of the invariant mass distribution after sub- 
tracting the combinatorial background. The raw spectra 
are corrected for the detector response (energy resolu- 
tion), the reconstruction efficiency, and occupancy effects 
(e.g. overlapping clusters). These corrections are made 
by embedding simulated single tt^'s from a full GEANT 
simulation of the PHENIX detector into real events, and 
analyzing the embedded tt" events with the same analysis 
cuts as used with real events. 

After computing the invariant yields in Au-|-Au col- 
lisions the medium effects are quantified using the 
nuclear modification factor (Raa)- Raa is the ratio 
between the measured yield and the expected yield for 
point-like processes scaled from the p-l-p result, and is 
defined as: 



IIaa(pt) = 



{l/Nlf)d^NfjdpTdy 
{TAA)d'^cr^N/dpTdy 



where (Taa) is the average Glauber nuclear overlap func- 
tion for the Au-f-Au centrality bin under consideration 



AA 



/TAA(b)db 



J(l_e-'^««^AA(b))^b 



/ AfO-jvivX /_mel 
W^coll lI'^NN 



(1) 



(2) 

where {N^^" ) is the average number of inelastic nucleon- 
nucleon collisions for the Au-|-Au centrality bin under 
consideration calculated with inelastic nucleon-nucleon 
cross section crJ^^J^- 

The measured nuclear modification factors for 0-5% 
central Au-|-Au reactions at ^/SNN = 200 GeV are shown 
in Fig. [T] and tabulated in Table |T] with statistical and 
systematic uncertainties. The systematic uncertainties 
fall into three categories. We denote Type A uncer- 
tainties as systematic uncertainties which are uncorre- 
lated from point-to-point. Type A systematic uncertain- 
ties are added in quadrature with statistical uncertainties 
and are shown as uncertainty bars. Partially correlated 
point-to-point systematic uncertainties are broken into a 
100% correlated component, referred to as Type B, and 
the above mentioned Type A. The Type B uncertain- 
ties are shown as gray boxes. The sources of type B 
uncertainties are discussed in detail in 9] , and are dom- 
inated by energy scale uncertainties, but also have con- 
tributions from photon shower merging at the highest pT 
{rv 15 — 20 GeV/c). There are also global systematic un- 
certainties, referred to as type C uncertainties, that are 
globally correlated systematic uncertainties (i.e. where 
all data points move by the same multiplicative factor). 
The type C uncertainties are ±12% and derive from un- 
certainties in the calculated nuclear thickness function 
and from the p-l-p absolute normalization. 



III. THEORETICAL CALCULATIONS 

The Parton Quenching Model (PQM) [3| encodes the 
dissipative properties of the system in terms of a single 
transport-coefficient, often referred to as q, obtained as 
the product of the parton-medium cross-section times the 
color-charge density. The average q quantifies the aver- 
age squared transverse momentum transferred from the 
medium to the parton per mean free path. The PQM 
model is a Monte Carlo program constructed using the 
quenching weights from BDMPS 0, [H. BDMPS is 
a perturbative calculation explicitly including only co- 
herent radiative energy-loss for the parton via gluon 
bremsstrahlung. The PQM model incorporates a real- 
istic transverse collision geometry, though with a static 
medium. It is also notable that the PQM model does 
not include initial state multiple scattering or modified 
nuclear parton distribution functions (PDF's). 

The Gyulassy-Levai-Vitev (GLV) @ model is a 
formalism developed to calculate in-medium gluon 
bremsstrahlung. An analytic expression is derived for 
the single gluon emission spectrum to all orders in opac- 
ity, assuming an infrared cut-off given by the plasmon 
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TABLE I: The tt^ nuclear suppression factor Raa as a function of transverse momentum for 0-5% Au+Au collisions at y/s^ n — 200 GeV. The 

type A, B, and C uncertainties are tabulated for each point. 



PT (GeV/c) 


Raa 


Type A Uncertainty 


Type B Uncertainty 


Type C Uncei 


1.25 


0.347 


±0.007 


±0.033 


±0.041 


1.75 


0.398 


±0.007 


±0.040 


±0.047 


2.25 


0.387 


±0.007 


±0.042 


±0.046 


2.75 


0.289 


±0.006 


±0.032 


±0.034 


3.25 


0.235 


±0.005 


±0.027 


±0.028 


3.75 


0.21 


±0.005 


±0.024 


±0.025 


4.25 


0.198 


±0.005 


±0.024 


±0.023 


4.75 


0.193 


±0.006 


±0.023 


±0.023 


5.25 


0.172 


±0.006 


±0.021 


±0.020 


5.75 


0.180 


±0.007 


±0.021 


±0.021 


6.25 


0.171 


±0.007 


±0.020 


±0.020 


6.75 


0.189 


±0.007 


±0.022 


±0.022 


7.25 


0.184 


±0.008 


±0.022 


±0.022 


7.75 


0.179 


±0.008 


±0.021 


±0.021 


8.25 


0.178 


±0.010 


±0.021 


±0.021 


8.75 


0.170 


±0.011 


±0.020 


±0.020 


9.25 


0.180 


±0.014 


±0.022 


±0.021 


9.75 


0.226 


±0.019 


±0.028 


±0.027 


11.00 


0.190 


±0.014 


±0.026 


±0.022 


13.00 


0.153 


±0.020 


±0.027 


±0.018 


15.00 


0.329 


±0.063 


±0.065 


±0.039 


17.00 


0.264 


±0.093 


±0.065 


±0.031 



frequency. Thus, within this framework, one can extract 
the local color-charge density. The color-charge density 
is written simply as dN^ /dy, assuming a completely glu- 
onic medium, or an equivalent dN'^'^ /dy for a mixture of 
quarks and gluons. In [12], using a realistic transverse 
collision geometry, the authors calculate a priori, with- 
out energy-loss, the single fixed geometrical average path 
length from the production point to the medium edge, 
{L)prod, and use it to calculate the parton energy-loss in 
a Bjorken expanding medium [l3| . The calculation also 
incorporates initial state multiple scattering effects and 
modified nuclear PDF's. 

The Wicks-Horowitz-Djordjevic-Gyulassy (WHDG) 
model utilizes the generalized GLV formalism [l^ for 
radiative energy-loss described above. In addition, their 
calculation includes a convolution of radiative energy-loss 
and coUisional energy-loss mechanisms. A realistic trans- 
verse collision geometry with a Bjorken time expansion 
is utilized, and then a full distribution of parton paths 
through the medium is calculated. The WHDG model 
does not yet include initial state multiple scattering or 
modified PDF's. 

The Zhang-Owens- Wang- Wang (ZOWW) [7] calcula- 
tions incorporate a next-to-leading order perturbative 
QCD parton model with modified jet fragmentation func- 
tions. The calculation explicitly includes only radiative 
energy-loss. A hard-sphere transverse collision geometry 



with a one dimensional expanding medium is utilized. 
The calculation also incorporates initial state multiple 
scattering effects and modified nuclear PDF's. 

The top left panel of Fig. [5] shows a comparison of 
the experimental data with calculated results from the 
PQM energy- loss model (as described in [4]) correspond- 
ing to different (q) values [IQ]. Note that only a subset 
of all the calculations corresponding to different (q) val- 
ues are shown in the figure for clarity. The upper 
right panel shows the suppression factor predicted at 
Pt = 20 GeV/c from the PQM model as a function of 
the (q) value. One can see that as the (q) increases, the 
additional suppression becomes smaller (i.e. saturates). 
This saturation effect was noted in |l3|, and interpret- 
ted as a result of the dominance of preferential surface 
emission. 

The other panels of Fig. [2] show similar comparisons 
of the experimental data with calculated results utilizing 
the GLV numerical calculation framework (as described 
in 12]) corresponding to different dN^/dy values [l7|; 
the WHDG calculational framework (as described in Q) 
corresponding to different dN^/dy values [13]; and the 
ZOWW calculational framework (as detailed in [tI]) cor- 
responding to different eo values [13] ■ Note that all calcu- 
lations are shown only for pt> 5 GeV/c as that is where 
the calculations are considered applicable. 

In Fig. [21 the same ir'^ suppression factor predicted 
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FIG. 2: (Color online) Left panels show it° Raa for 0-5% Au+Au collisions at ^^sWW=200 GeV and predictions from PQM 3, GLV O, 
WHDG [1, and ZOWW models with (from top to bottom) (q) values of 0.3, 0.9, 1.2, 1.5, 2.1, 2.9, 4.4, 5.9, 7.4, 10.3, 13.2, 17.7, 25.0, 40.5, 
101.4 GeV^/fm; dN<i /dy values of 600, 800, 900, 1050, 1175, 1300, 1400, 1500, 1800, 2100, 3000, 4000; dN" /dy values of 500, 800, 1100, 1400, 
1700, 2000, 2300, 2600, 2900, 3200, 3500, 3800; and cq values of 1.08, 1.28, 1.48, 1.68, 1.88, 2.08, 2.28, 2.68, 3.08 GeV/fm. Red hues indicate the 
best fit cases of (top) (q)= 13.2, (upper middle) dN^ /dy = 1400, (lower middle) dN" /dy = 1400, and (bottom) cq = 1-88 GeV/fm. Right panels 

show i?AA at pt — 20 GeV/c. 
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FIG. 3; (Color online) The nuelear suppression factors at pT — 20 
GeV/c for PQM as a function of (q)are shown as a blue line with a 
log-x and log-y display. Also shown is the functional form A 
(q)/(q> ~ 2.0 • ARaa/Raa) over the range 5 < (q) < 100 GoV^/fm, 



at pt= 20 GeV/c from the PQM model as a function 
of the (q) value is shown, but in this case with a log- 
X and log-y scale. The Raa ~ 0.75/\/(^ with (q) in 
units of GeV^/fm over the range 5 < (q) < 100. This 
means that over this range for a given fractional change 
in (q) there is always the same fractional change in Raa 
(i.e. A (q)/(q> « 2.0 • ARaa/Raa). 



IV. COMBINED RESULTS 

The task is now to detail how the experimental uncer- 
tainties (statistical and systematic) constrain the model 
parameters that are reflected in the suppression factors. 
As described previously, the uncertainties of the mea- 
sured points are separated into Type A (pT-uncorrelated, 
statistical ® systematic, Ci), Type B (px-correlated, abi, 
boxes on Fig. 1), Type C (normalization, uniform frac- 
tional shift for all points, CTc), where the cr's represent the 
standard deviations of the assumed Gaussian distributed 
uncertainties. With the predicted theory value fii{p) for 
each data point calculated for different values of the in- 
put parameter p, we perform a least-squares fit the to the 
theory by finding the values of p, ei,, tc that minimize: 



X {£b,ec,p) 



E 



iVi + eb'Jbi + ecUiCTc - p-i{p)f 



(3) 



where eb and tc are the fractions of the type B and C 
systematic uncertainties that all points are displaced to- 
gether and where ct^ = ai{yi + e^ai,^ + ^cViCrc)/Vi is the 
point-to-point random uncertainty scaled by the multi- 
plicative shift in Hi such that the fractional uncertainty is 
unchanged under systematic shifts, which is true for the 
present measurement. For clarity of presentation, the 
derivation of Eq. [3] (above) is given in Appendix A. 

For any fixed values of Cc, Eq. [3] follows the 
distribution with n + 2 degrees of freedom, for testing 
the theoretical predictions ^i{p), because it is the sum 
of n -t- 2 Gaussian distributed random variables. The 
best fit, Xmin' the minimum of x^(ef,,ec,p) by variation 
of ec, and p, is found by standard methods (for exam- 
ple using a MINUIT type minimization algorithm) and 
should follow the distribution with n ~ 1 degrees of 
freedom. The correlated uncertainties of the best fit pa- 
rameters are estimated in the Gaussian approximation 
by x^{eb,ec,p) = xLin + ^'^^ ^ standard deviation 
(cr) uncertainties. 

The present experimental type B uncertainties have 
point-to-point correlations whose exact correlation ma- 
trix is difficult to evaluate precisely. Thus, we consider 



two limiting correlation cases. The first is where the type 
B uncertainties are 100% correlated, i.e. all points move 
in the same direction by the same fraction of their re- 
spective type B uncertainty. The second is where the 
type B uncertainties are correlated such that the low px 
and high px points may shift with opposite sign (and 
linearly scaled in between), thus tilting the Raa either 
upward or downward as a function of pt . The minimum 
of the two cases is used for each constraint. 

We take as a first example the resulting theory predic- 
tions from the PQM model. For each calculation char- 
acterized by the (q), we calculate x^. We determine this 
value by varying et and Cc (the systematic offsets) until 
we obtain the minimum x^. These values are shown in 
Fig. m in the top panel. One can see that the overall 
lowest x^ value corresponds to (q) « 13 GeV^/fm. 

We then calculate the p- value for the x^ (the minimum 
of Eq. [3]), where the p- value is defined as: 



p- value = / Xin^){z)dz (4) 
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PQM Model ( q > (GeV^/fm) 



FIG. 4: (Color online) The statistical analysis results from the 
comparison of the PQM model with the tt^ Raa{pt) experimental 
data. The top panel shows the modified for different values of the 

PQM (q) . The middle panel shows the computed p-valuc directly 
from the modified as shown above. The bottom panel shows the 
number of standard deviations (cr) away from the minimum (best) 
(q) parameter value for the PQM model calculations. 




GLV Model dl^/dy 

FIG. 5: (Color online) The statistical analysis results from the 
comparison of the GLV model with the tt^ Raa{pt) experimental 
data. The top panel shows the modified for different values of the 
GLV dN^ /dy. The middle panel shows the computed p-value directly 
from the modified as shown above. The bottom panel shows the 
number of standard deviations {a) away from the minimum (best) 
dN^ /dy parameter value for the GLV model calculations. 



where xfna) chi-square distribution with the ap- 

propriate number of degrees of freedom, nj. This 
calculation is valid since the goodness-of-fit statistic 
follows a standard distribution. Note that p-value is 
the probability, under the assumption that the hypoth- 
esis is correct, of randomly obtaining data with a worse 
fit to the hypothesis than the experimental data under 
test [131 • These p- value's are shown in the middle panel 
of Fig. 12 

We find the overall minimum (or best) x^, and then 
compute which (q) scenarios are one and two standard 
deviations away from this minimum. The PQM trans- 
port coefficient (q) is constrained by the experimental 
data as 13.2 \ and ^ GeV^/fm at the one and two 
standard deviation levels, respectively. The two standard 
deviation constraints correspond to the 95% confidence 
intervals. We note that this range of large (q) values 
is currently under intense theoretical debate (see for ex- 
ample (U)- Thus, the quoted (q) constraint is for the 
model-dependent parameter of the specific PQM imple- 
mentation, and relating this parameter to the fundamen- 
tal value of the mean transverse momentum squared ex- 
change per unit length traversed may substantially lower 
the value. 



We apply the identical procedure to the GLV, WHDG, 
and ZOWW calculations and show those results in 
Figs. [HI [SI and [T] For the GLV calculations this results 
in a constraint of dN^/dy = lAQQt™ and ^290 at the 
one and two standard deviation levels, respectively. Con- 
straints for the WHDG model are dNs/dy = MOOtg^g 
and ^540 at the one and two standard deviation lev- 
els, respectively. Constraints for the ZOWW model are 
eo = 1.9lgJ and ^.^'^q at the one and two standard devi- 
ation levels, respectively. All of these constraint results 
are summarized in Table [Til 

For each of the above fits, there is a best fit value of 
£(, and Ec corresponding to the parameters in Eq. [3l For 
completeness, these values are for PQM, GLV, WHDG, 
and ZOWW, = 0.6 and = -0.3, ^ 0.7 and = 
—0.0, £(, = 2.1 and ec = —1.5, and ti, — 1.1 and tc = 
—0.6, respectively. All of the models considered here have 
a steeper px-dependence of Raa than the experimental 
data. Thus, the best fit is obtained within the type B 
uncertainties by tilting the Raa- 

It is notable that although there is a well defined over- 
all minimum in the modified for all four models, the 
maximum p-value in each case is different. In the PQM, 
GLV, and ZOWW models the maximum p- value's are ap- 
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WHDG Model dr^/dy 

FIG. 6: (Color online) The statistical analysis results from the 
comparison of the WHDG model with the tt^ i?AA(PT) experimental 
data. The top panel shows the modified for different values of the 
WliDG dN^ /dy. The middle panel shows the computed 
p-value directly from the modified as shown above. The bottom 
panel shows the number of standard deviations (ct) away from the 
minimum (best) dN^ /dy parameter value for the WHDG model 
calculations. 



proximately 9.0%, 5.5% and 7.8%, respectively. However, 
in the WHDG model the maximum p-value is substan- 
tially smaller at 1.3%. This is due to the fact that the 
WHDG model has a steeper px-dependence of the nu- 
clear modification factor Raa (regardless of parameter 
input) than the other models, and also steeper than the 
experimental data points. 

The identical best value of dN^ jdy for the GLV and 
WHDG calculations is interesting, since the inclusion 
of important coUisional energy-loss in WHDG leads to 
the naive expectation that a smaller color-charge den- 
sity would be needed for a similar suppression. How- 
ever, the different treatments of the distribution of paths 
through the medium may be compensating for this ef- 
fect. Further theoretical studies are needed to disentan- 
gle the physics implications of the shape differences and 
the similar best dN^/dy from GLV and WHDG. For aU 
the models considered, the relevant parameter constraint 
is approximately a ±20 — 25% uncertainty at the one 
standard deviation level. 

It is also interesting to inquire what simple linear fit 
function best describes the experimental data for pT> 5 
GeV/c. The identical procedure to that described above 




ZOWW Model (GeV/f m) 

FIG. 7: (Color online) The statistical analysis results from the 
comparison of the ZOWW model with the tt*^ /?aa(pt) experimental 
data. The top panel shows the modified x^ for different values of the 
ZOWW cq. The middle panel shows the computed p- value directly 
from the modified x^ as shown above. The bottom panel shows the 
number of standard deviations (tr) away from the minimum (best) eo 
parameter value for the ZOWW model calculations. 



is applied to the function Raa{pt) = b-|-m- px to deter- 
mine the best values for the two parameters. The best fit 
line and the envelope of lines with one standard deviation 
uncertainties are shown in Fig. [5] The results including 
all types of uncertainties are 6(intercept) = O.lGSlg 032 
and m(slope) = 0.0017t5]:E]E]|^ (c/GeV). The uncertain- 
ties on these parameters are correlated as shown by the 
one, two, and three standard deviation contours in Fig.[51 

Thus the data are consistent with a completely flat p^- 
dependence of i?AA for pt> 5 GeV/c (i.e. to = 0) within 
one standard deviation uncertainties. The maximum p- 
value for this simple linear function fit is 11.6%. 

The p-value's for all models considered are less than 
12%. It is notable that the five highest pt points 
(pt> 9.5 GeV/c), contribute over 70% to the total . 
As a check on the influence of these points on the ex- 
tracted parameter values, we have repeated the above 
procedure to the restricted range 5 < pt< 9.5 GeV/c. 
We find the following new constraints: PQM model (q) = 
13.2 t\% GeVVfm; GLV model dNS/dy = 1400l^?J]; 
WHDG model dN^/dy = IQQQt^l; ZOWW model 
eo = l-5lo 2 GeV/fm; simple linear fit ^(intercept) = 
0-170io:o34 and m(slope) = 0.00131°™^^ (c/GeV). We 
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TABLE II: Quantitative constraints on the model parameters from the PQM, GLV, WHDG, and ZOWW models and a linear functional form fit. 



Model 
Name 


Model 
Parameter 


One Standard Deviation 
Uncertainty 


Two Standard Deviation 
Uncertainty 


Maximum 
p- value 


PQM 


(a) = 13 2 GeV^ /fm 


+2.1 —3.2 


+6.3 —5.2 


9.0% 


GLV 


dNO /dy = 1400 


+270 -150 


+510 -290 


5.5% 


WHDG 


diVS /dy = 1400 


+200 -375 


+600 -540 


1.3 % 


ZOWW 


eo = 1.9 GeV/fm 


+0.2 -0.5 


+0.7 -0.6 


7.8 % 


Linear 


b (intercept) = 0.168 


+0.033 -0.032 


+0.065 -0.066 


11.6% 




m (slope) = 0.0017 (c/GeV) 


+0.0035 -0.0039 


+0.0070 -0.0076 
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FIG. 8: (Color online) The tt^ nuclear suppression factor i^AA as a 
function of transverse momentum for 0-5% Au+Au collisions at 
y^Siviv— 200 GeV. Point-to-point uncorrelated statistical and 
systematic uncertainties are shown as uncertainty bars. Correlated 
systematic uncertainties are shown as gray boxes around the data 
points. The global scale factor systematic uncertainty is quoted as 
text. Also shown are the best fit and the envelope of lines with one 
standard deviation uncertainty for a simple linear fit function 
constrained by the statistical and systematic uncertainties. 
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FIG. 9: (Color online) Shown arc the best fit values for m(slope) and 
b(intercept) as constrained by the experimental data. Also shown are 
the one, two and three standard deviation uncertainty contours. 



find that the resulting new constraints are within the one 
standard deviation uncertainties of those quoted for the 
full pt range. However, with the restricted range, the p- 
value's increase to 55%, 36%, 17%, 62%, and 75% for the 
PQM model, GLV model, WHDG model, ZOWW model, 
and the simple linear fit, respectively. Improvements in 
the data for pT > 9.5 GeV/c expected from future mea- 
surements will be crucial in determining whether any of 
the models discussed provide a statistically valid descrip- 
tion of the data over the full range 5 < pt < 20 GeV/c. 



v. SUMMARY AND CONCLUSIONS 

In this paper, we have compared model predictions 
of parton energy-loss with experimental data of semi- 



inclusive single high transverse momentum 7r° suppres- 
sion in central Au+Au reactions at y^SNN=200 GeV. In 
the comparison, statistical and systematic uncertainties 
were taken into account. We have obtained experimental 
constraints on model parameters of the color-charge den- 
sity of the medium or its transport coefficient. These 
values indicate a large medium density. It is crucial 
to note that the quoted constraints on these parame- 
ters do not include any systematic uncertainties in the 
models, but rather give the limits assuming a "perfect 
theory" with one unknown parameter, for example the 
color-charge density, constrained by the measurements 
including the experimental statistical and systematic un- 
certainties. Additional theoretical systematic uncertain- 
ties from the time evolution, energy-loss approximations, 
and calculation details need further investigation. 
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APPENDIX A: CONSTRAINT FORMALISM 

In the case of only point-to-point uncorrelated uncer- 
tainties (statistical and/or systematic), if one assumes 
they are Gaussian distributed and characterized by a, the 
root-mean square (RMS), calculating the best-parameter 
fit is straightforward via a log-likelihood or least squares- 
method 

The likelihood function C is defined as the a priori 
probability of a given outcome. Let 2/1,2/2, be n 

samples from a population with normalized probabil- 
ity density function f{y,p) where p represents a vec- 
tor of k parameters. For instance j/i could represent a 
measurement of a cross section at transverse momen- 
tum (pT)ij where the probability density of the measure- 
ment is Gaussian distributed about the expectation value 

= (y>: 



fiy,p} 



1 



cxp 



21 



(Al) 



If the samples are independent, then the likelihood func- 
tion is: 



C = l[f{y:,P) 



« exp- 



.1=1 



2<Tf 



(A2) 



However, if the samples are correlated, for example via correlated systematic uncertainties, then the full covariance 
matrix must be used 



Vij =< (y* - M»(p))(% - A*j(p)) > 

Then the likelihood function takes the more general form: 



C 



\V\ 



exp 



j^j^ {y,~^^,(p})V^-\y^~^^^{p)) 
1=1 j=l ^ 



(A3) 



(A4) 



where \V\ is the determinant of the covariance matrix V. Note that Eq. IA4I reduces to Eq. IA2I if the correlations 
vanish so that the covariances are zero and Vij is diagonal 

Vtj =< {vi - ^it{p}){yj - Mi(p)) % < (y, - fiiip)f >= Sijaf . (A5) 

Since Gaussian probability distributions are inevitable (as a consequence of the Central Limit Theorem) and since 
there is also an important theorem regarding likelihood ratios for composite hypotheses, it is convenient to use the 
logarithm of the likelihood 



2 In £ = In |U| + n log 27r + ^ ^ (y, - fiMW^j'iy, ~ H {p)) 

i=i j=i 



(A6) 
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We separate the uncertainties into four classes: type 
A) point-to-point uncorrelated systematic uncertainties; 
type B) correlated systematic uncertainties, for which 
the point-to-point correlation is 100% by construction, 
since the uncorrelated part has been separated out and 
included in uncertainty A); type C) overall systematic un- 
certainties by which all the points move by the same frac- 
tion (i.e. normalization uncertainties); and type D) sta- 
tistical. Categories A and D are simply added in quadra- 
ture and represent the total point-to-point uncorrelated 
uncertainties, denoted Ui below. 

We model a correlated systematic uncertainty as if 
there were an underlying uncertainty, e.g. absolute mo- 
mentum scale, which may cause correlated systematic 
variations Ay^ of the set of measurements, j/i, around 
their nominal value, that can be represented as a random 
variable, Z(,. The correlated type B variation of the mea- 
surements is represented by the displacement of all points 
from their nominal values by the correlated amounts 



{Ay''{sys).,Ay''{sys)j) 



(A7) 

;crfc,o-6., (A8) 



where Azf, = Zf, ~ {z^) — zi,. 

Since {zb) = 0, ((Azb)^) = (z^) _ (^,): 
random variable Zf, is the same for all i measurements 
while bi is a constant of proportionality which may be 



(T^ and the 



different for each i. We define = hab, where ±|crb. | is 
the systematic uncertainty bar shown on each point (gray 
box on each data point in Fig. [1] ) and where bi may be 
of cither sign, as it is possible that one point could move 
up while its neighbor moves down. The random variable 
Zb is assumed to have a Gaussian probability distribution 
f{zb), with r.m.s. ab 



1 



cxp - 



{Azbf 



2al 



(A9) 



The type C variation is independent of the type B 
variation. It is similarly assumed to be caused by an 
underlying random variable Zc that results in a system- 
atic displacement of the measurement by an amount Ay^ 
with 

Ay^isys^/y, = Az, (AID) 
{{Ay'^{sys\/y{){Ay'^{sys),/y,)) = ol (All) 

where by definition CTc is the same for all points. 

We then assume that the Likelihood function factor- 
izes as the product of independent Gaussian probabilities 
as in Eq. IA21 but that the distributions are correlated 
through their dependence on the random variables Zb and 



^ = ^J{yi,zb-,zc,v)!{zb)!{zc) 

i 

1 1 



CTia-2 • • • OnObOc ,/2ti 



<n+2) 



exp ■ 



^ {yi + bjAzb + yjAzc 

\i=l 



{Azbf , {Az,f 



2al 



2o-2 



(A12) 
(A13) 



To account for the type B systematic uncertainty, we allow any given sample of measurements, yi, corresponding to 
theoretical predictions fii (p) to have a correlated variation from their nominal values by an amount corresponding to 
a certain fraction eb of the underlying root-mean-square variation of zt, i.e. Azt = ebfJb, such that each point moves 
by an amount Ay^ — bietdt = ctUb^, the same fraction et of its systematic uncertainty bar; and similarly for the type 
C uncertainty. Then the likelihood function for any outcome, including the variation of and Sc would be: 



C = 



0-1(72 • ■ • (TnCTbCrc 



{n+2) 



exp ■ 



E 

^4=1 



2af 



2 



(A14) 



where the last two terms represent (Azb)2/(2a^) = e^cr^/(2cr^) and (Azc)^/(2cr^) = e^(T^/(2cr^) since we assumed the 
probability of the systematic displacements f{zb,c) to be Gaussian. Other probability distributions for the correlated 
systematic uncertainty could be used. For instance if ± jab; | had represented full extent systematic uncertainties, with 
equal probability for any Az;,, then the e^/2 term and associated normalization constant l/t7f,V27r would be absent 
from Eq. IA14I 

Then we use the likelihood ratio test to establish the validity or the confidence interval of the theoretical predictions 
/ii(p). One can use the modified log likelihood 



2 In/: 



■ icyiCTc - Mj(p))^ 



\i=l 



= X^(e6,ec,p) 



(A15) 
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because we will eventually take the ratio of the likeli- 
hood of a given set of parameters p to the maximum 
likelihood when all the parameters e;,, Ec and p are varied 
(the minimum value of Eq. lAlSp so that the terms pre- 
ceding the exponential in Eq. IA14I cancel because they 
are not varied. Eq. IA15I follows the x^-distribution with 
n + 2 degrees of freedom because it is the sum of n -I- 2 
independent Gaussian distributed random variables (i.e. 
in statistical terminology x^{eb, ^c,P) is xfn+2) )• This es- 
tablishes Eq. lAlSl as the x^-distributed quantity that we 
use for least squares fit to the theoretical predictions in- 
cluding the systematic uncertainties. Note that Eq. IA15I 
agrees with Eq. 8 in [23| in the discussion of fits with cor- 
related systematics. The specific procedure is described 
in the next paragraph. 

First 'fit the theory' to the data by minimizing Eq. lAlSI 
by varying all the parameters to find ic, p, the 
values of the parameters which give the overall mini- 
mum Xmin- If the Xmin ^'^^ ^his fit is acceptable for the 
n + 2 ~ (m + 2) = n — ni degrees of freedom, where 
TO are the number of parameters in then the the- 
ory is not rejected at this level. A confidence interval 
is then found for testing any other set of k parameters 
constrained to specific values, poj by again finding the 
minimum of Eq. IA15I for the k fixed values of po , by let- 
ting all the other parameters including Ch and Cc vary. For 
constant values of (7^, and large values of n, the "likeli- 



hood ratio" -2\Ti[C{pa)/C{p)] = -2[\nC{po) - \nC{p% 
i.e. x^(pb) ^ Xmin x^-distributed with k degrees of 
freedom, from which the confidence interval on the pa- 
rameters can be evaluated. However, in general, the un- 
certainty on the parameters is estimated in the Gaussian 
approximation by x^(Po) = Xmin + ^^ f'^'' ^ standard de- 
viation uncertainties (for example using a MINUIT (23j 
type fitting algorithm). 

For the present data, the statistical and random sys- 
tematic uncertainties are such that the shift in the mea- 
surement Hi due to the correlated systematic uncertain- 
ties preserves the fractional uncertainty. In this case the 
maximum likelihood and least squares methods no longer 
coincide and we use a least squares fit of Eq . [ATHl instead 
of Eg. [XTsl to estimate the best fit parameters: 




,(A16) 

where di is the uncertainty scaled by the multiplicative 
shift in Hi such that the fractional uncertainty is un- 
changed under shifts 



(Ji = O-'t 



(A17) 
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